Force network analysis of jammed solids 



S.-L.-Y. Xu\ X. Ilia 1 *, and J. M. Schwarz 1 
1 Physics Department, Syracuse University, Syracuse, NY 13244 
(Dated: August 27, 2010) 

Using a system of repulsive, soft particles as a model for a jammed solid, we analyze its force 
network as characterized by the magnitude of the contact force between two particles, the local 
contact angle subtended between three particles, and the local coordination number. In particular, 
we measure the local contact angle distribution as a function of the magnitude of the local contact 
force. We find the suppression of small contact angles for locally larger contact forces, suggesting the 
■ existence of chain-like correlations in the locally larger contact forces. We couple this information 

with a coordination number-spin state mapping to arrive at a Potts spin model with frustration and 
correlated disorder to draw a potential connection between jammed solids (no quenched disorder) 
and spin glasses (quenched disorder). We use this connection to measure chaos due to marginality in 
j^jy the jammed system. In addition, we present the replica solution of the one-dimensional, long-range 

, Potts glass as a potential toy building block for a jammed solid, where a sea of weakly interacting 

. spins provide for long-range interactions along a chain-like backbone of more strongly interacting 

^ ' spins. 
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I. INTRODUCTION 

The jamming phase diagram depicts a plausible scenario for a unified description of the phase change/crossover 
& ■ from a liquid to an amorphous solid in nonequilibrium systems ranging from granular particles to colloidal particles 
+^ ' to molecular particles These phase changes can be driven by temperature changes, packing fraction changes, 
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and/or changes in applied shear stress. Whether or not the boundaries of the jamming phase diagram are sharp 
in the equivalent equilibrium sense is still a matter of debate, particularly for molecular particles in the presence 
of temperature changes However, numerical studies of repulsive, soft particles at zero-temperature indicate the 
potential of a transition in the equivalent thermodynamic sense as the packing fraction of the system is increased 0, [f| . 

The numerics of the zero-temperature repulsive, soft particle system suggest that the transition from liquid to 
amorphous solid is somewhat of an unusual nature 043- For example, the average coordination number jumps from 
zero to some finite value at the transition, followed by a square root increase as a function of distance from the 
transition, i.e. <j) — <p c , where <f> denotes the packing fraction and <j) c , the critical packing fraction above which the 
^ ' system behaves as a solid. Moreover, the jump in the average coordination number only depends on the dimension 
of the system. This is reminiscent of the universality of the value of the jump in the spin-wave stiffness in the 
two-dimensional Kosterlitz-Thouless transition Q. While there exists a jump in a possible order parameter at the 
\ transition, the square root scaling demonstrates that the transition is not the typical first-order transition with no 
. ■ diverging correlation lengths. In fact, there are several diverging correlation lengths above and below the transition. 
For example, the fraction of jammed configurations as a function of <fi becomes increasingly sharper as the system size 
is increased. This sharpening can be tied to an increasing lengthscale below which the system acts as uncorrelated 
subsystems each with their respective critical packing fractions. 

As for the jammed phase itself, it exhibits various interesting properties. While the scaling of the bulk modulus 
with packing fraction is the same as ordered elastic solids, the scaling of the shear modulus with packing fraction 
shows an anomalous response 0, @] . In addition, the vibrational modes are of a "swirly" nature and appear to be 
quasi- localized as the jamming transition is approached from the solid side [tlHH- Similar behavior has been observed 
in other amorphous solids |ll| - [l4| . Both of these observations may be related to the non-affine displacements of the 
particles, which become more pronounced near the jamming transition. There is also an excess in the vibrational 
density of states beyond the Debye prediction at the transition. This excess has been likened to the Boson peak [l5l - fl7j 
in glassy systems. 

There have been a number of approaches to understand some of the above properties-jncluding a statistical mechan- 
ics approach [l8j . a field theoretic approach fl9| . and a phenomenological approach pol - |22j . The phenomenological 
approach examines the consequences of the isostatic nature of the transition by establishing a lengthscale below 
which the system behaves as an isostatic solid and above which the system behaves as an ordinary elastic solid. This 
lengthscale has been inferred from the numerical simulations by relating the frequency below which there exists an 
excess in the density of states (beyond the Debye prediction) to a lengthscale via a linear dispersion relation [23| . 
However, a direct measurement of the elasticity of the system via a local perturbation shows the system to respond 
on average as an ordinary elastic solid at all lengthscales at the jamming transition [24l. [25j. On the other hand, the 
study of fluctuations in the response of all contacts some radial distance from the local perturbation demonstrates a 
clear crossover lengthscale where the root-mean-squared fluctuations go from one type of power-law dependence on 
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radial distance to another type further away from the perturbation. This crossover lengthscale appears to scale with 
distance from the transition. 

One would like to understand why averaged response of the system behaves as ordinary elastic solid at the tran- 
sition. Or, conversely, why a potential isostatic lengthscale may be showing up only in fluctuations to a response 
in position space. In order to begin to investigate this mystery, we numerically study the force network of a two- 
dimensional system of repulsive, soft particles as the jamming transition is approached from above. In particular, we 
present evidence for the existence of chain-like structures in the locally large force bonds and argue that the particles 
surrounding the chain-like structures give rise to long-range interactions along the chain-like structures. However, 
the long-range interactions can be averaged out by the proliferation of these chain-like structures just as long-range 
interactions are typically screened out in polymeric systems. In addition, we propose a mapping of the coordination 
number to a spin state in order to potentially better understand some of the properties of the jammed solid from 
what will turn out to be a spin glass perspective. Therefore, we will be able to draw direct links between spin glasses 
and amorphous solids. Our mapping should further open the door for other amenable analytical techniques used in 
spin glasses with quenched disorder for studying amorphous solids where the disorder is not quenched. We note that 
there exists very interesting work on a replica-inspired approach to the glass transition for hard spheres as approached 
from the liquid side (as opposed to the solid side) [26j |. 

The paper is organized as follows. Section II presents our numerical analysis of the force network in the jammed 
phase. Section III establishes a connection of the jammed solid and spin glasses via analysis of the connectivity, or 
coordination number. Section IV goes beyond the numerical findings in Sections II and III to speculate on a plausible 
spin glass scenario for some properties of the jammed phase. Section V concludes the paper with a discussion of the 
relation of this work to preexisting work in the field. 



II. FORCE BONDS AND SPATIAL CORRELATIONS 



To study the force network of a two-dimensional jammed solid of repulsive, soft particles, jammed configurations are 
generated using the algorithm introduced by O'Hern and collaborators [E, @j- More specifically, the system consists 
of 50:50 mixture of N particles with a diameter, a, ratio of 1.4. The packing fraction sets the radii for a system of 
length unity. The particles interact via the following two-body potential: 



nry) = ~(i-^) A e(i-% (1) 

where is the distance between the centers of the two particles i and j, = (<7j +<7j-)/2, and e sets the energy scale 
for the system. The particles are placed randomly in the system with periodic boundary conditions and the conjugate 
gradient algorithm is invoked until the system reaches its nearest local minimum. We use an energy tolerance per 
particle of 10 -16 in units of e. 

In the repulsive, soft particle system, there are several known properties of the forces between overlapping particles: 
(1) the shape of the distribution of forces [27], [H| and (2) the distribution of forces demonstrates a lack of self- 
averaging [a, As for the shape of the distribution, the distribution is not long-tailed. Recent simulations in 
two- and three-dimensions have determined the shape of the tail down to normalized force values of 10~ 45 (28[. In 
two-dimensions, the tail is Gaussian. In three-dimensions, the tail falls off slightly more slowly than Gaussian, i.e. 
P(f) ~ eT^ ( j where / is the magnitude of the force between two overlapping particles. As for the lack of self- 
averaging, the distribution of forces takes on one form if the forces in each sample are normalized via the average 
force per sample or if the forces from all samples are pulled and normalized by the average force from all samples. 
Therefore, one cannot consider a large sample to be comprised of many smaller samples. A lack of self- averaging is 
also found in configurations of spin glasses. 

While the distribution of forces is a very useful quantity, it does not encode any spatial information. Motivated by 
the recent work of Zhou and Dinsmore [29|], we search for spatial correlations in the magnitude of forces. To do this, 
we measure the angle between any two force bonds emanating from a particle in the jammed packing at a particular 
packing fraction, <fi. We denote this angle as 9. See the schematic in Figure 1. The probability distribution for 8, 
P(0), is plotted in Figure 1 for all coordination numbers greater than 2 for <f> = 0.841 with A = 3/2 and N = 512. 
The lower bound on the coordination number is determined by the principle of local mechanical stability. For this 
particular <f> and system size, the average coordination of the jammed configurations, < z((f>) >= 4.076(2). Also, the 
fraction of jammed configurations is approximately two-thirds. So we are in the jammed phase. We should point out 
that earlier work extrapolated to a critical value of </> c in the infinite system limit of approximately 0.842. However, 
there exists a body of recent work suggesting that <p c depends on the protocol for obtaining jammed configurations 
such that "Point J" should be modified to "Segment J" (30l - [33| . Implications of these findings have not been fully 
fleshed out to date. Figure 1 also depicts P(0) for a particular z, or P z =i{ff). 
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FIG. 1: Left: P(0) and P(6i g ) for 4> = 0.841, N = 512, and A = 3/2. Right: P z =a{0) and P z = 4 (6»; 9 ) for the same parameters. 
The bin size is 2°. 



We then compare P(9) with P(6i g ), where 9i g is the angle between the two largest force bonds on a particle. We 
observe that the probability distribution in the latter case is more heavily biased towards the larger angles. In other 
words, there exists a suppression of the smaller angles between the two largest force bonds on a particle giving rise 

to chain-like correlations in the locally largest force bonds. Defining W — J 120 o (P(&lg) ~ P{9))d9 as a measure of the 
bias, W = 0.434(1) for = 0.841 and N = 512. 

Let us compare P{&) and P(9i g ) with the equivalent distributions for a more ordered packing. To do this, we use a 
monodisperse distribution of repulsive, soft particles in two-dimensions at a higher packing fraction. For a hexagonal 
packing of hard particles in two-dimensions, the packing fraction is equal to approximately 0.907. We choose a slightly 
larger packing fraction to ensure overlaps. See Figure 2. We observe three dominant peaks at 60°, 120°, and 180° 
degrees indicating a hexagonal packing. We also note a bias in P(9i g ) towards the larger angles for the ordered case 
as well. 

In the disordered case, there are a number of dominant peaks in P{ff) at = 54°, 60°, 64°, and 70°. There are 
also some subdominant peaks, many of them separated by intervals of 6°. See, for example, the range of 9 = 102° 
to 9 — 114°. These subdominant peaks are due to the bidispersity in particle radii. If the 1.4 ratio is increased, 
the spacings between the subdominant peaks increases. To see this, consider a particle surrounded by all "small" 
particles. If we replace one of the small particles by a larger particle, the larger particles will push its neighbors away 
due to insufficient space and therefore change the contact angle. We can calculate this change given the ratio between 
the two different radii in the just-touching case. For a radius ratio of 1.4, the change in contact angle is approximately 
5.4°, the interval between the subdominant peaks. 

For N — 512 and A = 3/2, as <f> is increased from cf> = 0.841 to cj) = 0.843, the bias towards the larger angles 
remains robust. More specifically, W remains constant. See Figure 3. For even larger (j>, W decreases. For example, 
W = 0.424(2) for = 0.846. This trend also agrees with the findings of Zhou and Dinsmore [29j where larger z 
suppresses chain formation. Since < z(cf>) > increases with 0, this trend is expected. Figure 4 also demonstrates this 
effect where P z =3A,5{^ig) arc plotted individually for <f> = 0.843. As z is increased from 4 to 5, W decreases from 
0.454(2) to 0.413(2). Finally, the trend of the suppression of smaller contact angles for larger force bonds persists 
for larger system sizes indicating that the suppression is not a finite system effect. We have also checked that the 
suppression persists for A = 5/2 and A = 2 as well, the more studied cases. 

To contrast our results with those of Zhou and Dinsmore [29j], we do not observe a peak in P z (0i g ) at 180 degrees 
for all z. Rather the peak in P z (9i g ) for z = 3, for example, is closer to 160°. Only for z = 4, is the peak at 180°. 
Therefore, the large force bond propagation is not as "straight" as in Zhou and Dinsmore [29[ simulations where a 
peak at 180° is observed for all z. Their protocol generates stable configurations by demanding force balance on each 
particle given that some of the forces on each particle are randomly generated. Zhou and Dinsmore argue that the 
finding of chain-like correlations is simply a consequence of Newton's third law. The peaks in their P z (9i g ) may be an 
indication of this. However, since our equivalent distributions are not all peaked at 180 degrees, a more complicated 
mechanism may be at work. 

Our results suggest a possible spatial correlation in the force bond strength. Given that larger angles between 
the locally largest forces persist, perhaps locally large force bonds can persist/percolate across the sample. To test 
for this, we construct a force chain by traversing along largest and second largest local force bonds and determining 
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FIG. 2: P(0) and P(6i g ) for (j> = 0.91, N = 512, and A = 3/2 for a monodisperse system. 



whether the chain of particles spans the system or not. The algorithm for this task is the following: 

(1) Identify the largest force bond in the system. 

(2) Traverse from one of the particles on either side of the largest force bond to the another overlapping particle 
along the second largest force bond. 

(3) Traverse to the next overlapping particle along either the largest force bond or the second largest force bond (if 
the largest force bond has already been traversed). 

(4) Repeat (3) until reaching an overlapping particle where both the largest and second largest force bonds have 
already been traversed. 

(5) Go back to the largest force bond in the system. Pick the other particle not initially chosen. If its second largest 
force bond has not been traversed, repeat steps (4) and (5). 

Note that the largest force bond of one particle is not necessarily the largest force bond of the other particle associated 
with the bond. This algorithm constructs the largest force chain in the sample. We also construct the weakest force 
chain in the sample by replacing largest with smallest, etc. Finally, we also study force chains beginning from any 
force bond in the system. 

In Figure 5 we present an example of the largest force chain. Loops make up some fraction of the "chain" . This is 
due to fluctuations in the forces. The larger force bonds exhibit chain-like structures, while smaller force bonds form 
local loops. We have checked this in the simulations by constructing the equivalent weakest force chain. Since the 
largest local force bond between two particles may be one of the smaller force bonds, as compared to the rest of the 
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FIG. 3: Left: P(9) and P(6 lg ) for <p = 0.843, N = 512, and A = 3/2. Right: P z =a{0) and P z=i (8 lg ) for the same parameters. 
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FIG. 4: Left: P z = 3 (0i a ) for 
for the same parameters. 



= 0.843, N = 512, and A = 3/2. Middle: P z = 4 (0i g ) for the same parameters. Right: P. 



force bonds in the system, chains can end in loops. Note that, at least for this example, the largest force chain spans 
the system in the vertical direction. 

How typical is this spanning for the largest force chains? Figure 6 shows the probability of spanning in either 
direction, P s , as a function of N for a particular 0. We see that P s decreases with increasing TV. While for small 
N, P s appears to decay linearly with N, the larger particle number data makes this candidate function less likely. 
We also observe that P s increases as <f> — > </> c , particularly for the largest particle number data, which is consistent 
with the increase in small angle suppression as <f> — > <j) c for the locally largest forces. As for the infinite system 
limit, it may be that the conditions to generate spanning chains need to be relaxed to generate spanning chains with 
nonzero probability. Also, if we start from any force bond in the system and generate a force chain moving along the 
locally largest two force bonds for each particle, the probability for spanning decreases. To be specific, for <f) = 0.841, 
N = 1024, and A = 3/2, P s w 0.563 starting at any force bond as compared to the largest force chain spanning of 
probability, where P s w 0.688. 

At this point, we must comment on the relation of our results to those of Makse and collaborators [34| performing 
dynamical simulations of deformable grains. They construct spanning force chains by starting with a grain at one 
end of the system and traversing the maximum force bond at every particle such that they reach the other side of 
the system. They find fewer spanning force chains closer to the transition than further away. They argue that this 
observation is due to an increasingly homogeneous distribution of forces far away from the transition. We should 
also contrast our results with work by Ostojic and collaborators (35[. They study the spatial extent of force bonds 
exceeding some threshold force. The threshold force is lowered until the force bonds exceeding the threshold force 
span the system. They observe a percolation transition of a new universality class. We, instead, investigate the local, 
largest forces. Also, we do not have spanning as a criterion, but rather as an "afterthought" . 

To better understand the spatial properties of the largest force chain, we measure its fractal dimension. Presumably, 
the fractal dimension is unity, though one should check this. To do so, we count the number of particles participating 
in the largest force chain. To relate this number with a length, in two-dimensions, L ~ y/N. Figure 7 plots the average 
number of particles participating in the largest force chain, Npc, as a function of system size, N . On the log- log scale, 
there is some curvature to the data. If we assume a fractal form, for A = 3/2 and <f> — 0.843, we measure a fractal 
dimension of 1.10(5)-very close to unity. Comparing our largest force chain data to the Ostojic formulation [35| . 
we measure the number of particles participating in the largest force spanning cluster, Npp. The fractal dimension 
of the spanning cluster of largest forces is 1.62(2). We then measure the number of particles participating in the 
smallest force spanning cluster, Ngp, where we replace "exceeding a threshold force" in the Ostojic formulation [35[ 
with "below a threshold force", the fractal dimension is 1.68(2). The fractal dimension of the spanning cluster at the 
ordinary percolation transition is about 1.89 [36j . which is somewhat larger than the value measured for the spanning 
cluster of largest forces (and smallest forces). This discrepancy may be due to the chain-like correlations in the largest, 
local force bonds, though one cannot discount the possibility of an eventual crossover to ordinary percolation. 



III. COORDINATION NUMBER AND SPINS 



The forces are not the only information encoded in the force network. The coordination number is another piece of 
information. Figure 8 is the same jammed configuration as depicted in Figure 5. Each color now represents a different 
coordination number. If the color of each vertex, a particle, shared by an edge, a contact, were of a different color, 
then this graph would represent a proper vertex coloring. However, note that there are some vertices sharing an edge 
that do have the same color. So Figure 5 does not represent a proper vertex coloring. 



FIG. 5: Jammed configuration for <f> — 0.841, N = 1024. The particles participating along the largest force chain (dark blue) 
as distinguished from the other particles (light blue). 
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FIG. 6: The probability of spanning in either direction, P s , for the largest force chain. 



But let us, for the moment, address proper vertex colorings. A proper vertex coloring using at most p colors exists 
if there exists a zero energy solution of zero-temperature antiferromagnetic Potts model with p states. In other words, 
the p-coloring problem maps to a zero-temperature antiferromagnetic Potts model with p states. Now, if one were 
to consider the jamming problem, the number of possible colors participating in the jammed clusters is 5 (z = 3 
through z = 7, for this particular two-dimensional bidisperse system) which leads to a p = 5 state Potts model. At 
<f> = 4> c , the isostatic condition would bias one of the colors — for example, the color associated with coordination 
number four dominates for two-dimensional frictionless discs. This bias can be encoded via a magnetic field. The 
isostatic condition also imposes other weights on the other colors as well such that the average coordination number 
is four. In particular, recent work shows that at the jamming transition in two-dimensions, the coordination numbers 
range from 3 to 5 such that z — 4 particles make up about half the system and z = 3, 5 particles make up the other 
half |37|- As the packing fraction increases, the weights for each color changes. For a hexagonal packing, all couplings 
are ferromagnetic. 

Both the antiferromagetism and the competing magnetic fields (to account for the average coordination number) 
contribute to frustration of the spin system. Moreover, there exists randomness in the system due to (1) the fact that 
not all particles are participating in the jammed configuration and (2) the randomness of the forces. The combination 
of frustration and disorder should lead to a spin glass phase [H, |39j]-a spin glass phase with no quenched disorder. In 
fact, the couplings and dilutions (particles not participating in the jammed configurations) are dynamically generated. 
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FIG. 7: Log-log plot of the average number of particles in the largest force chain that span, Nfc, the average number of 
particles in the largest force spanning cluster, Nlf, and the average number of particles in the smallest force spanning cluster, 
Nsf, as a function of N. The error bars are smaller than the symbols. The straight lines are fits to the data. Here, (f> — 0.843 
and A = 3/2. 



However, the lack of self-averaging in the forces is characteristic of quenched disordered system. Indeed, there are 
systems with no disorder that can be mapped to quenched disorder and there are systems with unquenched disorder 
that can also be mapped to quenched disorder. For the former, consider an Ising spin system whose interaction 
fluctuates in sign with distance between sites poj . For the latter, consider a liquid of hard spheres as described by 
hypernetted chain equations [26| . 

The potential for a jamming-spin glass correspondence provides motivation to analyze the jamming system in 
terms of a spin glass. If we compute the ratio of antiferromagnetic couplings (different coordination numbers (colors) 
between contacts) to ferromagnetic couplings (same coordination numbers (colors) between contacts), for cf> = 0.841, 
N = 1024, for example, the ratio is approximately 2.5. We could ignore the random variation in the magnitude of the 
forces (interactions between spins) and set the coupling to be a constant using the measured ratio of antiferromagnetic 
spins to ferromagnetic spins to simulate a two-dimensional p = 5 Potts glass system at zero-temperature. While there 
exists a number of ground state algorithms for two-dimensional Ising spin glasses (4lTj43| , we do not know of an exact 
ground state algorithm for the 5-state Potts glass and so we leave this avenue for future work. 

Instead, to further the possible connection between jammed solids and spin glasses, we compute the spin glass 
equivalent of bond chaos in the repulsive, soft particle system [iij . Chaos in mean- field spin glasses is presumably due 
to the fact that the spin glass phase is a marginal phase such that a perturbation in the energy of the system is sufficient 
to alter the weights of different equilibrium configurations [isf . More specifically, in mean-field, the Hessian associated 
with Parisi's ansatz for the structure of the matrix order parameter contains all non- negative eigenvalues (4fj| . In 
finite-dimensions, the marginality is presumably due to the slowly decaying correlation functions |47j |. It is interesting 
to note that jammed solids at the jamming transition are marginally rigid and, therefore, may also exhibit similar 
chaotic features. 

So, we use this notion of chaos in spin glasses to look for a quantitative measure of chaos in the jamming system. 
To do this, we first assign a list of random numbers as the initial position of particles in a system and use conjugate 
gradient algorithm to generate a jammed state. Next, we use the same initial positions and perturb them by a 
given strength, then apply the same conjugate gradient algorithm to generate another jammed state. Specifically, for 
particle i with initial position (xi, yi) in the unperturbed system, the corresponding particle i in the perturbed system 
has an initial position of 



where S is the magnitude of perturbation strength, and Ci and di are randomly generated numbers chosen from the 
same distribution used to generate Xi and yi. 

If the two sets of initial positions both lead to a jammed state, we calculate the overlap between the two states as 
defined by 



y'i = Vi + S d i 



(2) 
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FIG. 8: Jammed configuration (same as above) colored via coordination number. Light blue denotes z = 0, magenta denotes 
z — 3, red denotes z = 4, blue depicts z = 5, orange depicts z — 6, and purple denotes z — 7 (possible for 1.4 diameter ratio 
bidisperse system). 
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FIG. 9: Plot of r as a function of N for several 8s. Inset: Plot of r for N = 128 as a function of S. The error bars are smaller 
than the symbols. 



where Sj is the spin state of the unperturbed system (using the Domb representation), S'j is the spin state of the 
perturbed system, and M is the number of spins in both systems with z > 2, which grows with N. The brackets 
denote ensemble averaging. For S = 0, r = 1. The system is chaotic when 

lim lim r(S,N,(f>) < 1. (4) 

S— >0 N— >oo 

It may be that the limit of <j> — > 4> c may also have to be taken to observe chaos given that only at the jamming 
transition is the system marginally rigid. 

In Figure 9 we present results for rasa function of N for different values of 8 and <fi. For fixed <f>, we observe that r 
initially decreases with increasing TV and then begins to level off (though one cannot rule out an small increase in r for 
the largest N studied). Moreover, r increases with decreasing S for fixed N, as expected. Since it is not clear that r 
reaches a limit as a function of N for fixed 8, we cannot easily extrapolate to the zero perturbation limit. However, it 
is clear from the data that as (f> decreases towards unjamming, r decreases, indicating that as the jamming transition 
is approached from above, the system is becoming increasingly chaotic. 
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IV. JAMMED SOLIDS AND SPIN GLASSES 

Let us summarize our numerical findings: 

(0) There is a suppression of the smaller angles subtending locally larger force bonds. This suppression in- 
creases as 4> c is approached from above. 

(1) The probability for spanning of the largest force chain increases as <p c is approached from above. 

(2) The coordination number-spin state mapping suggests that the jammed solid becomes increasingly chaotic as <f) c 
is approached from above. 

Keeping the above findings in mind, we recall one of the traditional models of an amorphous solid via a randomly 
dilute network of springs with probability p, the model otherwise known as rigidity percolation [Hj]. While the 
nature of the rigidity percolation transition remains contentious [49, _5Q] , an important concept has emerged in terms 
of constraint counting, a concept initially presented by Maxwell back in 1864 [51] . Below the rigidity percolation 
threshold, there are only underconstraincd bonds, above the threshold, there exists overconstrained bonds. At the 
transition, on average, the system is not overconstrained or underconstrained — it is isostatic. This condition should 
then determine the location of the transition. For central force interactions, p r — 2d/z, where d is the dimension of 
the system and p r is the occupation probability above which the system is rigid. For the triangular lattice, p r — 2/3. 
Numerical studies find a result close to this estimate. However, the Maxwell argument does not take into account 
redundant bonds and the possibility of a fractal set of overconstrained bonds participating in a rigid cluster. 

Jacobs and Thorpe have numerically extended Maxwell's argument by identifying those overconstrained bonds 
via an algorithm implemented on a two-dimensional bar-joint network [52j, which makes uses of a theorem from 
Laman [53j j . They demonstrate that at the transition, the overconstrained bars make up a fractal subset of the bonds. 
As more bars are added, the fraction of overconstrained bars increases such that they make up a finite fraction of 
the bars. It turns out that the fraction of overconstrained bars can be viewed as the order parameter for the rigidity 
percolation transition. 

Indeed, it would be interesting to extend some of these ideas of rigidity percolation more concretely to the repulsive, 
soft particle system, such as trying to identify the overconstrained bonds via use of the pebble game algorithm. 
However, use of the pebble algorithm takes into account only fixed connectivity information and it does not take into 
account the force information. Of course, there have been several recent works drawing a more intimate connection 
between rigidity percolation and jamming. Ellenbroek and collaborators map the jammed solid to a network of 
stretched springs [54|. Study of the stretched spring network demonstrates anomalous scaling of the bulk modulus 
as opposed to the shear modulus. Even more recent work investigates a square lattice occupied by springs [55| . 
Such a lattice is isostatic with a frequency-independent density of states. As next-nearest springs are occupied, the 
system undergoes a rigidity percolation transition for an infinitesimal occupation of next-nearest neighbor springs. 
This system exhibits the kind of transition proposed by the phenomenological framework put forth by Wyart and 
collaborators. There exists a clear lengthscale below which the system behaves isostatically and above which it does 
not, which is not so readily apparent in the disordered packing, at least in terms of a direct measurement. 

While working within the more traditional spring framework is certainly useful, coupling the coordination number- 
spin state mapping with the force information via interaction strengths between the spins and the dilution information 
provides one with an alternate description of a jammed solid — a spin description to be compared to the traditional 
framework. For example, overconstrained bonds in rigidity percolation may correspond to frozen spins — spins that 
take on the same state in every configuration — in the spin picture. From the numerical information obtained so far, 
there are more antiferromagnetic interactions than ferromagnetic interactions. As for the magnitude of the interactions 
between the spins, recall the chain-like correlations observed in the locally larger force bonds. This information can 
be also be encoded into a Potts spin Hamiltonian via chain-like correlations in the spin interaction strengths. 

How does the Potts spin Hamiltonian change with increasing packing fraction? The amount of dilution decreases. 
Moreover, the ratio of antiferromagnetic to ferromagnetic interactions also decreases. What about the force informa- 
tion? We observe that the probability of spanning increases with decreasing 4> for the largest system sizes. Based on 
this observation, the chain-like correlations in the spin interaction strengths decrease as the packing fraction increases. 
All of these properties can be encoded into a Potts spin Hamiltonian. 

Given the interplay of frustration and randomness, we conjecture that a Potts spin analog of the jammed solid should 
exhibit spin glass behavior. A possible spin glass/jammed solid correspondence may help to explain the ordinary 
elasticity findings of Ellenbroek and collaborators at all lengthscales when measuring the system's averaged response 
to inflation of a central particle [24], [25[. For both the paramagnet and the spin glass, the average magnetization is 
zero, i.e. one cannot distinguish between the two phases. To observe the spin glass phase, one must measure higher 
order quantities such the Edwards- Anderson order parameter. Ellenbroek and collaborators did observe a deviation 
from ordinary elasticity when measuring fluctuations in the forces averaged over all contacts within a distance from the 
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perturbation [2J, |25[ . The lengthscale beyond which the fluctuations level off may correspond to a finite-dimensional 
spin glass phase, which exhibits very different properties than a mean-field spin glass. The replica, or mean-field, 
scenario [Ea . gives evidence for infinite number of pure thermodynamic states, while the low-dimensional, or 
droplet scenario [58|,[5{|, at least in the Ising case, argues for two groundstates, the same ones occurring in the Ising 
ferromagnet — all spins up or all spins down. In the droplet picture, an applied magnetic field is a relevant perturbation 
driving the system to a ferromagnetic phase, while in the replica scenario, an infinitesimally applied magnetic field 
does not destroy the spin glass phase. How such a lengthscale separating mean-field and low-dimensional behavior 
could come about, we propose in the following subsection. 

A. Simple example: Long-range, one-dimensional Potts glass 

While one would ultimately like to simulate a zero-temperature, two-dimensional Potts-glass with chain-like corre- 
lations in the couplings, for the time being we make use of a simple example, namely, the one-dimensional, long-range 
interacting p-state Potts glass with Gaussian distributed quenched bond randomness. There exists previous analysis of 
the one-dimensional, long-range interacting Ising spin glass. In light of our coordination number-spin state mapping, 
we extend these results to the Potts glass case. 

Why study a model with long-range interactions? The chain-like correlations in the forces may call for an effective 
theory centered on these chains with the non- force chain bonds providing for long-range forces along the chain. In fact, 
the existence of these chain-like correlations is reminiscent of a picture proposed earlier by Cates and collaborators [(5(| 
where the jammed solid is comprised of linear force chains and a sea of spectator particles modeled as an incompressible 
solvent. Instead of considering the force chain as a linear object supporting loads only along its own axis, consider the 
force chain as a polymer — a polymer embedded in solution, i.e. the other particles. Just as the hydrodynamics of the 
solution couple one part of the polymer with a distant part, resulting in long-range interactions along the polymer, 
the force chain would also experience long-range interactions [6lj . We note that the polymer- in-solution analogy has 
its shortcomings. A more accurate analogy might be a polymer embedded in another polymer-type medium. Also, 
resumably the particles participating in the largest force chain would no longer be part of the largest force chain upon 
shearing or other perturbations. In other words, there is particle conversion between the force chain and non-force 
chain particles. 

How would this effective force chain framework depend on the packing fraction? At the onset of jamming, presum- 
ably there exists an isotropic, fractal network of larger force spanning chains. As the packing fraction is increased, 
the larger force chains become shorter such that eventually there is no distinction between the chains and the other 
particles since the definition of a force chain would have to be relaxed in order to achieve spanning. So the lengthscale 
of the chain lattice shrinks as packing fraction is increased and the long-ranged interactions become screened out 
just in the case of many polymers in solution. At a particular packing fraction; while for lengthscales smaller than 
"chain" lattice, the system behaves long-ranged, or mean-field-like; for lengthscales larger than "chain" lattice, the 
long-range interactions are screened out and the system behaves as a low-dimensional system. Therefore, within 
one system, one can interpolate between mean-field and finite-dimensions due to the correlated heterogeneity in the 
forces/interactions. 

While a precise formulation of an effective force chain theory is still lacking, we discuss the possible ramifications 
of long-range interactions in a conventional spin glass system. We begin with the p-state Potts glass Hamiltonian in 
rf-dimensions, 

j.. p- 1 

i,j |Z - 71 a=l 

where Si is a p-state Potts spin in the simplex representation [(32| at site i and a denotes the interaction range. The 
distribution of the coupling constant is given by 

The replica trick is applied to the Hamiltonian to compute the free energy averaged over disorder by considering n 
identical replicas of the original system. The averaged free energy, F, is then given by 

~Z™ - 1 

F=-k B T\im , (6) 

n^>0 n 
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where 

p-i 



Z» = J([[ P(J«)A7«)Tr {a?} exp J £ £ fc^fep E 1 " 

After integrating over the disorder and applying the Hubbard-Stratonovich transformation, 

^ = / n d ®f ab ex p { - J v«&9&} Tr m ex p {fia^osa} > ( ? ) 



where the matrix i^T is defined by 



Kij ~ (k B T)*\i-j\*° (8) 



and ^ < with J = kgT g such that we set the ferromagnetic order parameter to zero. 

To perform the trace over the spins, we recall that the simplex representation of the p-state Potts model requires 



s 

Y,<4 = J*W-1, (9) 

a 

E< - °> 

s 

where in the replicated partition function each Sj n is represented by one of the e Q 's. We also define v a i,c — ^2 S e a e i> e c 
for convenience 63[. We expand the second exponential term in Eq. 7 to order Q 3 , and after computing the trace 
and re-exponentiation, we obtain 



Tr {s?} exp {s?XQf ab } - \ E (Qta b ) 2 + 1 E Q% b QUcQ 



■ya. 
i.ca 



1 £ + 0((QSJ)*). (10) 

After taking the continuum limit, we rewrite Z n in momentum space as 

Z* = J II dQ a J(q)cM-J d d qH{Q a J(q)}) (11) 



with 



F = 7 /E( r + ^"W(e)^f(-?)-^/ E (2T) d %l + 32 + 33)Q:f(9i)Q^(©)Q2a(93) 

J «o# "'919293 a ^/3/ 7 



-— u 
2 



"19293 ^ 



where r, u, and w are the coupling constants and J q = J We have omitted a q 2 Q^{q)Q^{—q) term. 

Performing the momentum shell renormalization group with one- loop corrections using the usual notation [fjij , we 
obtain 

r' = ( 2 b- d (r-(-2(p-l)v 2 + l8 P A (p-2) 2 u 2 )C(b)) 
v' = ( 3 b- 2d (v + (-3p + 4)v 3 D(b)) 
u' = ( 3 b- 2d (u + 36u 3 p\ P -'3) 2 D(b)) 

£ _ b l+d/2-r,/2 

(13) 
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where C(b) = f* /b (2ff)d(g( d 2o . g - d)+r)2 and D(b) = fy b j^yqm^ij- Since the renormalization group does not generate 
new long-range terms, t] = 2 + d — 2a to all orders. Therefore, for a < both cubic couplings are irrelevant and 
the critical behaviour is mean-field (65[ with Gaussian fluctuations. For a > the cubic terms become relevant and 
different critical behaviour is observed. 

The modified value of rj = 2 + d — 2a is a general feature of long-range interactions. The v and u terms are irrelevant 
at a < 2/3 in the one-dimensional case. Thus, the system behaves mean-field-like when 1/2 < a < 2/3, where the 
lower bound a < 1/2 assures the convergence of free energy. This result is a simple generalization of the previous 
p = 2 case (6(| E3] and the short-range case [68] . Modifying the range of the interaction corresponds to modifying 
the effective spatial dimension of the system from low-dimensional to mean-field. This property has been numerically 
verified in the p = 2 case and the 3-spin spin glass [H, Iz3|- in the jammed solid, the interaction range is presumably 
tuned by the presence of other spanning force chains (Potts chains). On lengthscalcs smaller than the chain lattice, 
the long-range interaction dominates and above this lengthscale, the short-range interaction dominates. 



V. DISCUSSION 

The existence of force chains have been demonstrated both experimentally [7lj and numerically [34{ . We numerically 
demonstrate the existence of chain-like correlations in the locally large forces in the repulsive, soft particle system. 
This finding may call for an updated effective force chain description with the non-force chain particles providing 
for long-range interactions along the force chains as the basis of the description. One can test for these long-range 
interactions by perturbing a force chain locally and observing the response further along the force chain. 

Experiments by Majumdar and Behringer 72] have shown an ordering of force chains in the presence of shear. 
In the absence of shear, the network of locally large force chains is isotropic. In the presence of shear, the force 
chains align along the shear direction. A reordering of the force chain network structure via non-affine displacements 
such that the system is inherently anisotropic may account for the anomalous shear modulus exponent found in the 
repulsive, soft particle system. One can test for this reordering in the repulsive, soft particle system by checking for 
spanning of the locally large forces in the direction along and perpendicular to the shear separately. 

We also propose a coordination number-spin state mapping. With this mapping, we find that the repulsive, soft 
particle system becomes increasingly chaotic as the the jamming transition is approached from above. More testing 
of this mapping is needed. For example, one can investigate the response of the system to shear and see how the 
coordination number changes throughout the system. A changing coordination number corresponds to a changing 
spin state. It may be that the framework of random spin wave theory (or some alteration thereof) can be used to 
understand the response of the jammed solid to various perturbations. 

Based on the fact that jamming and glass transition share some commonalities, via transitivity, one suspects that 
the spin glass transition should be related to the jamming transition. A particle (graph) coloring with non-trivial 
biasing to enforce the isostaticity condition at the onset of jamming may provide us with the third missing link. 
Whether or not the equivalent spin glass is mean-field or low-dimensional remains to be seen. Of course, both the 
marginality and the massless, or replicon, modes found in mean- field spin glasses make for a tantalizing match [731173]. 
A link between jammed solids and spin glasses also provides a connection to constraint satisfaction problems. This 
connection has recently been explored by Mailman and Chakraborty [75[ and also by Krazkala and Kurchan [7(1 [73] ■ 
Constraint satisfaction problems are at the heart of attempting to understand what makes problems solvable. 

While the phenomenological isostatic framework provides insight into the jamming transition, its success hinges on 
the measured exponents for system sizes ranging up to approximately 10000 particles currently. History teaches us 
that the reliance on numerical results for support may ultimately reveal oversimplified assumptions [78j . Hence, the 
importance of constructing well-defined, calculable toy models of jamming transition cannot be underestimated. Such 
is the case for fc-core percolation, the simplest model of rigidity percolation [79} . If one allows the fc-core percolation 
transition to take place as soon as a non-zero spanning fc-core cluster is allowed, i.e. minimal rigidity, then the mean- 
field exponents are presumably in the same category as those measured for the finite-dimensional repulsive, soft sphere 
system. Note that minimal rigidity here does not translate to an average bond occupation of four per site. While 
the fc-core constraint is indeed a local one, it contains the nonlocal feature of rigidity percolation-that the removal 
of a single bond can trigger the removal of an entire occupied cluster. Other spin systems, such as the dilute 3-state 
Potts antiferromagnet on the triangular lattice, also exhibit a similar nonlocal property [S0|] . There exists other work 
towards calculable models of jammed solids that should be pursued further. For example, see Refs. [8ll483j . 

In closing, given the connection between spin systems and repulsive, soft particle systems, we expect that construc- 
tion and study of new spin models with constraints inspired by jamming will pave the way for further insights into 
jammed solids. Should one of these candidate spin models exhibit spin glass properties, we may come full circle to a 
percolation view of jamming as has been developed for spin glasses 
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